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SUMMARY 

^»This  paper  uses  a  very  general  version  of  the  Iterative  Proportional  Fitting  Procedure  to 
develop  an  algorithm  for  estimation  in  simultaneous  logit  models.  The  algorithm  can  be  used 
for  any  loglinear  model  which  can  be  cast  in  the  form  of  simultaneous  logit  equations.  The 
principal  advantage  of  this  method  is  that  it  is  not  necessary  to  fit  parameters  associated  with 
the  sampling  constraints  and  thus  very  large  problems  can  be  attacked.  A  numerical  example 
using  GLIM  and  a  sample  GLIM  macro  are  included. 


KEY  WORDS:  Categorical  Data.  Loglinear  Models.  Iterative  Proportional  Fitting  Procedure. 
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1.  INTRODUCTION 

Simultaneous  logit  models  have  been  oft-proposed  but  have  suffered  from  a  lack  of  a 
suitable  and  easy  estimation  technique.  We  show  how  an  application  of  some  simple  ideas 
which  underly  the  Iterative  Proportional  Fitting  Procedure  (IPFP)  leads  to  a  simple  and 
worthwhile  algorithm.  The  same  approach  can  be  used  for  any  loglinear  model  and  is 
potentially  useful  for  large  contingency  table  problems 

Methods  for  analyzing  data  in  the  form  of  a  binomial  response  variable  with  one  or  more 
covariates  are  well  known.  A  common  method,  and  the  one  which  corresponds  to  a  loglinear 
model,  is  logit  regression.  There  are  numerous  programs  (e.g.,  BMDPLR  and  GLIM)  which 
will  find  maximum  likelihood  estimates  for  the  logit  regression  model.  When  the  response 
variable  is  multinomial,  rather  than  binomial,  the  situation  is  rather  different.  The  models 
which  correspond  to  logit  regression  posit  a  system  of  simultaneous  logit  models  (see.  e.g., 
Fienberg  (1980).  Fienberg  and  Mason  (1978)  and  Nerlove  and  Press  (1973)).  One  method  of 
finding  M.L.E.'s  in  this  case  is  to  transform  the  simultaneous  logit  model  into  an  equivalent 
loglinear  model.  This  approach,  however,  introduces  potentially  many  nuisance  parameters 
corresponding  to  the  "interactions"  among  the  independent  variables.  It  is  possible  that  the 
introduction  of  these  extra  parameters  will  make  the  model  too  large  to  be  fit  conveniently. 
To  overcome  these  problems  Fienberg  (1980)  suggests  fitting  continuation  ratios  instead  of  the 

> 

maximum  likelihood  analysis.  We  will  show  that  this  procedure  is  actually  the  first  step  in  an 
IPFP  and  that  the  full  IPFP  is  itself  relatively  easy  to  calculate.  We  also  give  a  simple  GLIM 
macro  for  this  analysis. 

■  The  method  is  developed  in  several  steps.  First  we  consider  a  multinomial  response  variable 

'  with  categorical  explanatory  variables.  We  then  generalize  this  to  continuous  explanatory 

•  variables  and  finally  consider  models  more  complicated  than  those  normally  discussed. 

i 

In  order  to  develop  the  method  we  need  a  brief  discussion  of  some  theory  underlying  the 
IPFP. 


! 


3 


2.  BACKGROUND  AND  NOTATION 

Csiszar  (1976)  presents  a  very  elegant  and  general  discussion  of  the  IPFP  by  developing  a 
"geometry"  for  the  information  measure.  A  simplified  version  of  the  chief  results  of  this 
theory  are  outlined  below.  Let  n,  p.  q,  r,  s,  and  t  denote  p.m.f.'s  which  are  non-zero  fOT  all 
elements  of  a  finite  set  I.  The  Kullback-Leibler  information  number  (or  directed  divergence) 
specifies  a  distance. 

I(p||q)  =  yyy  I.  e  ,  P(0  In  (p(i)/q(i)) 

between  p  and  q.  The  principle  of  minimum  discriminant  information,  as  formulated  by 
Kullback  (1959),  aims  to  minimize  the  distance  between  a  reference  distribution,  q  above,  and  a 
family  of  other  distributions.  The  properties  of  such  estimates  have  been  studied  extensively. 
The  most  important  results  can  be  found  in  Kullback  (1959)  and  are  summarized,  with  a 
special  emphasis  on  contingency  tables,  in  Gokhale  and  Kullback  (1978).  Darroch  and  Ratcliff 
(1972)  also  used  the  directed  divergence  in  their  development  of  the  Generalized  Iterative 
Scaling  algorithm. 


We  need  to  specify  an  appropriate  family,  E,  of  probability  mass  functions  (p.m.f.'s)  over 
which  to  minimize  the  distance.  Linear  sets  of  p.m.f.’s  are  a  natural  class  and  correspond  to 
usual  loglinear  models.  A  convex  set,  E,  of  p.m.f.’s  is  called  linear  if  when  p  and  q  are  in  E 
and  t  =  a  •  p  +  (1-c)  •  q  (a  e  P)  is  a  p.m.f..  then  t  is  also  in  E.  A  p.m.f.  which  satisfies 


min  1  (p||r) 
p  e  £ 


is  called  the  I-projection  of  r  on  E  and  will  be  denoted  by  q  =  P^ir).  Csiszar  gives 
conditions  under  which  P[(r)  exists  (it  is  always  unique)  and  develops  a  geometry  for  I- 
projections  by  using  an  analogue  of  Pythagorous’  Theorem.  Now  let  F  =  (f  :  y  c  D  be  a  set 
of  real  valued  functions  on  I  and  A  =  (a^:  y  e  D  be  real  constants.  Define  Mp  to  be  span 
(F).  A  linear  set,  E,  can  be  constructed  by  considering  the  set  of  p  for  which. 


Z  p(i)  •  Mi)  =  a.  ;  y 
i  c  i  y 


e  r  . 
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When  we  consider  s  to  be  an  observed  probability  function  and 

a  =  I  s(i)  f  (i)  :  y  c  r 
y  i  e  1  r 

then  the  duality  between  maximum  likelihood  and  minimum  discriminant  (see  e.g.  Darroch  and 
Ratcliff  (1972))  estimation  states  that  if 

q  =  Pe(T) 

then 

ln(q)  e  Mp  +  ln(r) 
and 

x 

q  -  s  e  Mp  , 

i.e.  q  is  the  m.l.e.  (under  Poisson  sampling)  for  the  corresponding  log-affine  model.  Csiszar’s 
principal  theorem  says  that  if  £  is  the  finite  intersection  of  the  linear  sets  E  (i.e.  £  =U  E  ) 

k  c  k 

Le  K 

then  q  =  Pg( r)  is  the  pointwise  limit  of  qn  =  Pg  (q  )  n  =  1.2.3  .  .  .  where  q  =  r  and  En  = 

n 

£  _  if  i  =  n  mod  j  K  (  . 

It  is  this  theorem  on  cyclic  projections  that  we  shall  use  for  the  simultaneous  logit  algorithm. 
For  more  details  about  the  above  theory  see  Csiszar  (1976)  or  Meyer  (1981a)  and  for  another 
application  of  the  theory  see  Meyer  (1981b).  We  now  present  a  short  example  to  illustrate  the 
notation. 

Example  1 


Let  p  be  an  observed  3x3  probability  function  obtained  via  multinomial  sampling  and 
consider  the  model 

£(p(ij»  =  q(i,j>  , 

and  ln(q(i,j))  =  p  +  a.  +  p.  ;  ij  =  1,2,3,  i.e.  independence  of  row  and  column  categories. 


The  linear  manifold  for  this  model  is  spanned  by  a  set  of  tables,  f'  and  f1  ,  where 

1  k  =  i 


»•»  *  { 


0  k  *  i 


and 


r  (k.E)  =  { 


1  i  =  j 
0  £  #  j 


Now  lei 


F,  1  - 

Fc  *  {  f: :  j  *  1-2-3  } 


F  -  F.Uf,  . 


We  further  allow 


-  Ki  p<k-«  f^-«- 

aj.=  2a  p<k.f>  f^(k.O 


and  define  Er  and  Ec  to  be  the  linear  spaces  generated  by  the  fR.  aR  and  fc>  ac  pairs. 

The  maximum  likelihood  estimate  q(i,j)  can  be  generated  by  taking  q  (i,j)  =  1  Vi,j  and 
successively  forming  the  I-projections  onto  Er  and  £c.  This  algorithm  is  just  the  usual  IPFP 
and  converges  after  one  step. 

3.  SIMPLE  SIMULTANEOUS  LOGIT  MODELS  FOR  A  THREE-WAY  TABLE 
In  this  section  we  consider  simultaneous  logit  models  for  three  way  tables  where  there  is  one 
response  variable  and  two  explanatory  variables.  In  the  next  section  we  consider  the  general 
simultaneous  logit  case. 

Consider  a  multinomial  response  variable  A  (with  1  levels)  and  explanatory  variables  B  (with 

J  levels)  and  C  (with  k  levels).  Denote  the  table  of  observed  counts  by 

w  (i  =  1,...,  I;  j  =  1.....  J;  k  =  1 .  K) 

ijk 

where  w  is  the  number  of  responses  at  level  i  when  the  explanatory  variables,  B  and  C.  are 
at  levels  j  and  k  respectively.  In  this  discussion  we  transform  the  data  w  into  a  probability 
vector,  z..k  =  w  /  w  This  is  not  necessary  in  practice.  We  will  denote  the  expected  cell 
probability  £( z..k)  by  mjjk  . 

When  I  =  2  the  saturated  logit  regression  model,  using  the  notation  of  Fienberg  (1980).  is 
logit  (m  )  =  v  +  v  +  v  +  v 

6  ijk  2(j)  3<k)  23(jk) 

If  we  leave  out  the  interaction  term,  v2J,  the  model  becomes 
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logit  (mijkl  =  v  ♦  v,(ji  ♦  vj(t) 

with  likelihood  equations 


m  =  z 
i  r  >r 


“Vk  = 


and  the  sampling  constraint 

m  =  z  . 

’ji¬ 
lt  is  a  simple  matter  to  rewrite  this  model  as  a  loglinear  model,  viz.,  ln(m  , )  =  u  +  u  +  u, 

r  ijk  l(i>  2(  j 1 


+  u  +u  +u  +  u  with  likelihood  equations 

3(k>  Uliii  13<ikl  23<jk) 


m  =  z 

ij-  ij* 


m  =  z 

-jk  -jk 

The  logit  model  was  defined  only  for  I  =  2  but  the  loglinear  model  is  defined  for  any  value 
of  1.  If  we  consider  this  arbitrary  loglinear  model  we  can  write  it  as  the  set  of  simultaneous 
logit  models 

(3.1-D  V  =  v' +  i,/  VL 


(3.1-1) 

(3.1-2) 


ln(iVV =  y:  +  V  iik, 


(3.1— (1—1)) 
(3.1— (I)) 


= yl" + ij  * i;:, 

i“<m  /m  >  =  v1  ♦  v;(*  4, 


The  last  equation.  (3.1— (I))  is  redundant  but  it  doesn’t  hurt  to  include  it  and  it  may  actually 
improve  the  resulting  algorithm.  We  shall  similarly  break  up  the  likelihood  equations  into 
in,  =  z,  ,  in,  =  z, 

ij*  ij*  2r  2j+ 


(3.2-1) 


"Vk  =  Z,.k  ’  "Vk  =  Vk 


"Vk  +  K  =  Z.jk  +  Z2jk 
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.  K  =  Z:r  ' 

K  =  Z3.r 

(3.2-2) 

Ki  =  Z:-k 

•  ®3*k  =  Za-k 

K  +  : 

=  z  +  z 

3jk  3jk 

Ki.r  =  z.-i 

.  m,  =  z, 

•j*  ij*  ij* 

(3.2— (I— 1)) 

m,  =  z 

l-l,*k  1-1 

_3> 

♦ 

r~ 

It 

7- 

m,  , .  +  m, 

1-1 .jk  Ijk 

=  z,  +  z, 

1-1,  jk  Ijk 

K  =  V  • 

m  =  z 

ij-  ij- 

(3.2-(I)) 

Ki =  Vk  • 

m  =  z 

i-k  i-k 

m,  +  m  = 

Ijk  Ijk 

=  2  +  Z 

Ijk  tjk 

To  fit  the  full 

loglinear  model 

would  require  IJK  - 

very  large  number. 

Each  of  the 

logit  models  requires 

w 


needed  to  fit  only  the  logit  models  then  the  problem  would  be  easy. 


Consider  the  following  algorithm.  Let  m°k  =  1,  V  ijk.  Estimate  inland  mousing  the  logit 
model  (3.1-1)  with  likelihood  equations  (3.2-1)  and  set  m1  =  m°  for  i  *  1.2.  Now  estimate 

ijk  ijk 

m^  and  m;  from  (3.1-2)  and  (3.2-2)  with  m:  =  in'  for  i  *  2,3.  Continue  cvcling  through  the 
logit  equations,  each  time  adjusting  one  pair  of  the  multinomials  leaving  the  rest  untouched. 
Each  iteration  of  the  algorithm  is  a  (small)  logit  regression  problem  with  an  entire  cycle 
requiring  1  logit  regressions.  We  will  now  show  that  in'-  converges  to  the  M.L.E.  m  . . 

ijk  ijk 


If  we  view  the  problem  from  the  I-divergence  point  of  view  we  know  that  we  require 

m  =  P^(  m°) 

where  E  =  (p.m.f.’s  p ,  :  p  =  z. 

l  *  ijk  rij*  ij+ 

P  =  Z 

ri*k  i+k 

*.*  *  z..«' } 

If  we  view  the  logit  regressions  as  merely  calculating  I-projections  we  note  that 
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m1  =  f*£  (m°) 

u  i 


where 


and 


where 


etc. 


E(  =  {p.m.f.’s  p  satisfying  (3.2-1)} 
m:  =  (m1) 

E,  =  (p.m.f.’s  p^  satisfying  (3.2-2)) 


As  O'  E  =  E  (in  fact  D1'1  £  =  £)  we  can  invoke  the  theorems  on  I-projections  and  state 

i*l  1  1-1  l 

~  £ 

that  m  must  converge  to  the  M.L.E. 

l.lk 

Thus  the  procedure  of  sequentially  solving  the  logit  equations  is  really  a  version  of  the  1PFP. 
There  is  nothing  special  about  the  particular  logit  formulation  we  have  used.  All  that  is 
required  is  that  the  E  corresponding  to  the  logit  models  must  have  an  intersection  which  is  E. 
We  note  that  some  of  the  orderings  suggested  by  Fienberg  (1980.  p.  110)  do  not  satisfy  this 
requirement,  although  they  are  sensible  for  other  reasons. 


Generally  it  is  not  possible  to  obtain  estimates  of  asymptotic  (co)variance  matrices  using  an 
IPFP  algorithm.  In  this  case,  however,  if  the  method  used  to  solve  the  individual  logit 
problem  will  yield  asymptotic  (co)variance  estimates  then  this  information  can  be  used  to 
compute  (co)variance  estimates  for  the  multinomial  logit  problem.  Fot  example,  if  estimates 

for  the  v*"  and  their  covariance  matrices  for  l  =  1 .  I  are  obtained  on  the  last  cycle  of  the 

algorithm,  these  will  be  correct  for  the  full  problem.  It  is,  however,  impossible  to  find 
estimates  of  the  correlations  between,  say.  v|  and  v*,  for  exactly  the  same  reasons  that 
covariances  are  unavailable  in  the  usual  IPFP. 


In  order  to  clarify  the  above  discussion  we  present  a  small  numerical  example  of  a  trinomial 
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logit  regression  problem. 
Example  2 


Consider  two  explanatory  variables,  sex  and  intelligence,  each  with  two  levels.  For  each 

combination  of  levels  we  observe  a  trinomial  response.  The  data  are 

Response 
D  E  F 


1(1) 

S(l) 

7 

7 

n 

i 

S(2) 

8 

24 

20 

1(2) 

S(l) 

9 

6 

16 

S(2) 

26 

11 

16 

The  simultaneous  logit  model  we  consider  is 

ln(D/E)  -  G1  -  P  *  S' 

*  i 1  «.!' 

ln(E/F)  =  G:  -  1:  -  S: 

1 1  1 1* 

and  ln(F/D)  =  G  *  V  *  S? 

U'  l.|! 

Each  of  these  logit  models  requires  3  parameters  whereas  the  equivalent  loghneaT  model  has 
10  parameters  (6  corresponding  to  two  logit  models  and  4  from  the  sampling  constraints).  By 
using  a  loglinear  model  algorithm  we  obtain  the  M.L.E.  for  the  cell  means  as  approximately 

Response 


D 

E 

F 

Kl) 

S(l) 

4.035 

7.634 

9.331 

S(2) 

10.965 

23.366 

17.668 

1(2) 

S(l) 

11.965 

5.366 

13.668 

S(2) 

23.035 

11.634 

18.331 

After  one  cycle  of 

the  simultaneous 

logit  approach  (i.e.. 

3  logit  models)  the  fitted  values  are 

Response 

D 

E 

F 

1(1) 

S(l) 

4.217 

7.163 

9.621 

S(2) 

10.783 

23.837 

17.379 

1(2) 

S(l) 

11.783 

5.837 

13.379 

S(2) 

23.217 

11.163 

18.621 

Trjr » 


I 
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and  arte-  '  *  - . .c\  nt  •  uted  values  are  ihe  same  as  those  from  the  loglinear  model  algorithm. 

D  p 

If  ^ e  constrain  tne  parameters  so  that  r  and  S^.  are  zero,  for  l  =  1.2.3  .  then  the 
estimated  parameter  salues  and  their  standard  errors  are: 


Parameters 

Estimated  parameter 

Estimated  s.e. 

G 

-.638 

.4706 

[;, 

-.119 

.4737 

s' 

V 

1.440 

.4327 

'  G: 

-.201 

.4015 

4 

.480 

.4278 

4 

-.734 

.4023 

G1 

.839 

.4341 

4 

-.361 

.4111 

4 

-.705 

.4074 

While  it  is  possible  to  calculate  the  estimated  covariance  between,  say.  G1  and  from  the 
Newton  algorithm  used  for  the  logit  regressions,  it  is  not  possible  to  estimate  the  covariance 
between,  say,  G1  and  G:. 

The  calculations  for  this  example  were  performed  using  GLIM  and  the  complete  output  and 
macros  are  in  the  Appendix. 

4.  GENERAL  SIMULTANEOUS  LOGIT  MODELS 
Thus  far  we  have  only  touched  the  surface  of  models  which  can  be  fit  using  this  technique. 

but  the  general  form  is  the  same  as  the  specific  example  we  have  considered.  Let  J  =  11.2 . 

J}  be  an  index  set  and  suppose  that  for  each  element  of  J  we  observe  an  1-nomial  response. 
Denote  the  response  by  x  .  i  =  1 I,  J  e  J  and  the  expected  count  by  m  .  In  other  words 

'J  ij 

we  have  product-multinomial  sampling  with  x  observations  in  category  j.  Any  of  the  usual 
factorial  loglinear  models  which  condition  on  the  explanatory  variables  can  be  written  in  the 


11 


simultaneous  logit  form, 

ln(m)  /m, )  c  M 
lnlm,  /ir>3 )  c  M 

lnin^  /m^  e  M 

where  M  is  a  linear  manifold  in  R3.  The  linear  manifold  M  determines  a  linear  space  of 
p.m.f.’s  E.  Just  as  in  the  previous  case,  cycling  through  the  logit  equations  leads  to  the 
M.L.E.’s  for  the  corresponding  loglinear  model. 

In  the  above  formulation  there  is  no  need  for  each  logit  to  belong  to  the  same  linear 
manifold.  In  other  words,  the  most  general  form  of  loglinear  model  which  satisfies  the 
sampling  constraints  can  be  written  as 

lnin^  /m, )  f 
ln(m,  /m3 )  c  M, 

ln(m  /m  )  e  M 

M  l.  l-l 

where  each  M  C/7J.  Each  M  determines  a  linear  space  of  P.D.'s,  E  a  subset  of  the  p.m.f.'s  on 
I  x  J  points,  and  sequentially  fitting  the  logit  models  produces  the  1-projection  onto 
E  =  (>;  E .  Thus  any  loglinear  model  which  involves  product-multinomial  constraints,  or  can 
be  written  in  this  form,  can  be  fit  as  a  sequence  of  logit  models.  In  this  way  it  is  possible  to 
fit  the  loglinear  model  without  necessarily  fitting  "parameters"  for  the  product-multinomial 


constraints. 
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APPENDIX  I.  A  GLIM  MACRO  AND  AN  EXAMPLE 


1.1.  GLIM  MACROS 

The  following  GLIM  macros  were  used  to  calculate  the  multinomial  logit  estimates  of 
Section  3. 

$SUB  MLOG 

$PRI  '* DELETE  logithelp  IF  NOT  NEEDED  *  ’ 

SMAC  LOGITHELP  ! 

There  are  6  macros.  3  which  include  parameter  estimates  and  three  which  don’t.! 

WOUT  -  Basic  macro  with  parameter  estimates,  requires  macros  DM  and  PARA.! 

NOUT  -  Basic  macro  without  estimates,  requires  macros  DD  and  LOG.! 

I 

Typically  one  would  use  NOUT  until  satisfied  that  the  process  has  converged! 
and  then  WOUT  ,  once,  to  get  parameter  estimates.! 

NOUT  and  WOUT  have  up  to  9  arguments,  one  for  each  response  level.! 

Initialize  the  n  with  $ARG  NOUT  A  B  ...  I  .  where  A  B  ...  I  are  identifiers! 
for  the  response  levels.  The  scalar  %N  must  be  set  to  the  number  of  ! 
response  levels  and  NOUT  and  WOUT  must  have  %N  arguments.! 

To  get  started,  one  needs  to  fit  the  model  to  some  data,  any  will  do.! 

The  identifiers  A  B  ...  I  are  overwritten  with  the  fitted  values.  ! 
so  copv  them  or  vou'll  loose  them.! 

SENDM 


SMAC  NOUT  ! 

SCAL  %D=0  ! 

SOUT  SARG  LOG  %1  %2  %3  %4  %5  %6  %7  %8  %9  SUSE  LOG  ! 

SUSE  DD  SOUT  5  SPRI  ’CONVERGENCE  ’  %D  SEND 
SMAC  DD  ! 

SYVAR  %1  SCAL  7Z=%\*%2  SERR  B  ZZ  SFIT  .  SCAL  %D=%IF(%GT(%DV.%D).%DV.%D)! 
SCAL  %1=%FV  SCAL  %2=ZZ-%1  SEND 
SMAC  LOG  ! 

SARG  DD  %1  %2  SUSE  DD  SARG  DD  %2  %1  SCAL  %M=%EQ(%N.2)  SEXIT  %M  ! 

SARG  DD  %2  %3  SUSE  DD  SARG  DD  %3  %1  SCAL  %M=%EQ(%N.3)  SEXIT  %M  ! 

SARG  DD  %3  %4  SUSE  DD  SARG  DD  %4  %1  SCAL  %M=%EQ(%N,4)  SEXIT  %M  ! 

SARG  DD  %4  %5  SUSE  DD  SARG  DD  %5  %\  SCAL  %M=%EQ(%N.5)  SEXIT  %M  ! 

SARG  DD  %5  %6  SUSE  DD  SARG  DD  %6  %1  SCAL  %M=%EQ(%N.6)  SEXIT  %M  ! 

SARG  DD  %6  %7  SUSE  DD  SARG  DD  %7  %1  SCAL  %M=%EQ(%N.7)  SEXIT  %M  ! 

SARG  DD  %7  %8  SUSE  DD  SARG  DD  %8  %1  SCAL  %M=%EQ(%N.8)  SEXIT  %M  ! 

SARG  DD  %8  %9  SUSE  DD  SARG  DD  %9  %\  SEND 

SMAC  WOUT  ! 

SCAL  %D=0  ! 

SARG  PARA  %1  %2  %3  %4  %5  %6  %7  %8  %9  SUSE  PARA  ! 

SUSE  DM  SPRI  ’CONVERGENCE  ’  %D  SEND 
SC  need  to  get  routines  for  goodness  of  fit  here. 

SMAC  PARA  ! 

SARG  DM  %1  %2  SUSE  DM  SARG  DM  %2  %1  SCAL  %M=%EQ(%N.2)  SEXIT  %M  ! 
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$ARG  DM  %2  %3  $USE  DM  $ARG  DM  %3  %1  SCAL  %M=%EQ(%N.3)  SEXIT  %M  ! 

$ARG  DM  %3  %4  $USE  DM  $ARG  DM  %4  %1  $CAL  %M=%EQ(%N.4)  SEXIT  %M  .' 

SARG  DM  %4  %S  SUSE  DM  $ARG  DM  %5  %1  SCAL  %M=%EQ(%N,5)  SEXIT  %M  ! 

SARG  DM  %5  %6  SUSE  DM  SARG  DM  %6  %1  SCAL  %M=%EQ(%N.6)  SEXIT  %M  ! 

SARG  DM  %6  %7  SUSE  DM  SARG  DM  %7  %1  SCAL  9bM=%EQ(%N,7)  SEXIT  %M  ! 

SARG  DM  %1  %8  SUSE  DM  SARG  DM  %8  %1  SCAL  %M=%EQ(%N,8)  SEXIT  %M  ! 

SARG  DM  %8  %9  SUSE  DM  SARG  DM  %9  %1  SEND 

SMAC  DM  ! 

SYVAR  %1  SCAL  ZZ=%l+%2  SLRR  B  ZZ  SFIT  .  SCAL  %D=%IF(%GT(%DV,%D).%DV<%D)! 

SDIS  E  SCAL  %1=%FV  SCAL  %2=ZZ-%1  SEND 

SRETURN 

SFINISH 


1.2.  AN  EXAMPLE  USING  THE  MACROS 
The  complete  GLIM  output  for  Example  2  is  presented  below. 

GLIM  3.09  (01977  ROYAL  STATISTICAL  SOCIETY.  LONDON 

7SINPUT  35  MLOG  SC  ...  INPUT  THE  MACROS  $ 

♦DELETE  logithelp  IF  NOT  NEEDED  » 

7SPRINT  LOGITHELP  $ 


There  are  6  macros.  3  which  include  parameter  estimates  and  three  which  don’t. 

WOUT  -  Basic  macro  with  parameter  estimates,  requires  macros  DM  and  PARA. 
NOUT  -  Basic  macro  without  estimates,  requires  macros  DD  and  LOG. 

Typically  one  would  use  NOUT  until  satisfied  that  the  process  has  converged 
and  then  WOUT  .  once,  to  get  parameter  estimates. 

NOUT  and  WOUT  have  up  to  9  arguments,  one  for  each  response  level. 

Initialize  the  run  with  SARG  NOUT  A  B  ...  I  ,  where  A  B  ...  1  are  identifiers 
for  the  response  levels.  The  scalar  %N  must  be  set  to  the  number  of 
response  levels  and  NOUT  and  WOUT  must  have  %N  arguments. 

To  get  started,  one  needs  to  fit  the  model  to  some  data,  any  will  do. 

The  identifiers  A  B  ...  I  are  overwritten  with  the  fitted  values, 
so  copy  them  or  you’ll  loose  them. 

7SDEL  LOGITHELP  SUNIT  4  SDATA  D  E  F  SREAD 

1  1  1  1  8  24  20  9  6  16  26  11  16  $ 

?$F ACTOR  I  2  S  2  SCAL  I=%GL(2,2)  :S=%GL(2.1)  $ 

7SCAL  Tl=l  :T2=2  SYVAR  Tl  SERR  B  T2  SC  ...  GET  THE  ALGORITHM  STARTED  $ 
7SFIT  S*I  S 


SCALED 


CYCLE  DEVIANCE 
4  0.8195E-15 


1 


DF 


?$CAL  A=D  :B=E  :C=F  SC  ...  COPY  THE  DATA  SO  THAT  IT  ISN’T  LOST  $ 

?$CAL  %N=3  SARG  NOUT  ABC  SUSE  NOUT  $ 

CONVERGENCE  2.791 
7SLOOK  ABCS 


1 

4.217 

7.163 

9.621 

2 

10.78 

23.84 

17.38 

3 

11.78 

5.837 

13.38 

4 

23.22 

11.16 

18.62 

?SC  ...  THESE  ARE  THE  ESTIMATES  AFTER  THE  FIRST  ITERATION  $ 

?$USE  NOUT  $ 

CONVERGENCE  0.0950 

7$  USE  NOUT  $ 

CONVERGENCE  0.0013 

?$USE  NOUT  $ 

CONVERGENCE  0.0000 

? 

?$C  ...  NOW  USE  WOUT  TO  GET  PARAMETER  ESTIMATES.  ETC.  $ 
7SARG  WOUT  ABC  SUSE  WOUT  $ 

SCALED 

CYCLE  DEVIANCE  DF 

3  0.2651E-06  1 

ESTIMATE  S.E.  PARAMETER 

1  -0.6377  0.4705  %GM 

2  -0.1189  0.4737  S<2) 

3  1.440  0.4327  1(2) 

SCALE  PARAMETER  TAKEN  AS  1.000 

- CURRENT  DISPLAY  INHIBITED 

SCALED 

CYCLE  DEVIANCE  DF 

3  0.8541E-07  1 

ESTIMATE  S.E.  PARAMETER 

1  -0.2008  0.4015  %GM 

2  0.4803  0.4278  S<2) 

3  -0.7342  0.4023  1(2) 

SCALE  PARAMETER  TAKEN  AS  1.000 
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- CURRENT  DISPLAY  INHIBITED 

SCALED 

CYCLE  DEVIANCE  DF 


3  0.3044E-07 

1 

ESTIMATE 

S.E. 

PARAMETER 

1  0.8385 

0.4341 

%GM 

2  -0.3614 

0.4111 

S(2) 

3  -0.7054 

0.4074 

1(2) 

SCALE  PARAMETER  TAKEN  AS  1.000 

- CURRENT  DISPLAY 

INHIBITED 

CONVERGENCE 

0.0000 

?$LOOK  ABCS 

1  4.035 

7.634 

9.332 

2  10.97 

23.37 

17.67 

3  11.97 

5.366 

13.67 

4  23.03 

11.63 

18.33 

7 


?$C  ...  ONE  COULD  NOW  WRITE  MACROS  TO  TEST  GOODNESS-OF-FIT  OF 
ABC  WITH  D  E  F  S 

?$STOP 
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